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Abstract 

Here we present a detailed analysis of solar acoustic mode frequencies and their 
rotational splittings for modes with degree up to 900. They were obtained by applying 
spherical harmonic decomposition to full-disk solar images observed by the Michelson 
Doppler Imager onboard the Solar and Heliospheric Observatory spacecraft. Global 
helioseismology analysis of high-degree modes is complicated by the fact that the 
individual modes cannot be isolated, which has limited so far the use of high-degree 
data for structure inversion of the near-surface layers (r > 0.97-R©). In this work, we 
took great care to recover the actual mode characteristics using a physically moti- 
vated model which included a complete leakage matrix. We included in our analysis 
the following instrumental characteristics: the correct instantaneous image scale, 
the radial and non-radial image distortions, the effective position angle of the solar 
rotation axis and a correction to the Carrington elements. We also present variations 
of the mode frequencies caused by the solar activity cycle. We have analyzed seven 
observational periods from 1999 to 2005 and correlated their frequency shift with 
four different solar indices. The frequency shift scaled by the relative mode inertia 
is a function of frequency alone and follows a simple power law, where the exponent 
obtained for the p modes is twice the value obtained for the / modes. The different 
solar indices present the same result. 

Keywords: Helioseismology, Observations; Instrumental Effects; Oscillations, So- 
lar; Solar Cycle, Observations 



1. Introduction 

The central frequencies of solar acoustic modes, which are obtained using spherical 
harmonic decomposition, have been successfully used to determine the solar interior 
structure to as close as 21 Mm to the solar surface (r < 0.97i?©) using modes 
with angular degrees £ < 300 (e.g., Gough, Kosovichev, Toomre, et at, 1996). The 
inclusion of high-degree modes (i.e., up to £ = 1000) has the potential to improve 
dramatically the inference of the sound speed and the adiabatic exponent (ri) in 
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the outermost 2 to 3% of the solar radius, allowing to construct localized kernels as 
close to the solar surface as 1.75 Mm (Rabello-Soares, Basu, Christensen-Dalsgaard, 
et al., 2000). The effects of the equation of state, through the ionization of hydrogen 
and helium, are felt most strongly in the outer layers of the Sun, making this shallow 
region of particular interest. Furthermore, dynamical effects of convection, and the 
processes that excite and damp the solar oscillations, are predominantly concentrated 
in this region. Although the spatial resolution of the modern helioseismic instruments 
allows us to observe oscillation modes up to £ = 1000 and higher, only a small 
fraction of them are currently used (£ < 300). Unfortunately, analysis of high-degree 
data is complicated by the fact that the individual modes cannot be isolated (e.g., 
Rabello-Soares, Korzennik and Schou, 2001). 

The solar structure is not static, but changes over the solar cycle. It is well known 
that the mode frequencies change with solar activity. It seems that the responsible 
mechanism is restricted to the outer layers of the Sun (Libbrecht and Woodard, 
1990), where the high-degree modes are confined. However, at the moment, there 
is no general agreement in the precise physical mechanism that gives rise to the 
frequency variation. It is likely a product of the change in the subphotospheric small- 
scale magnetic field strength with the solar activity cycle (e.g., Goldreich, Murray, 
Willette, et al., 1991). Accordingly to Dziembowski and Goode (2004), the frequency 
shift is easily explained in terms of a variation in the turbulent velocities associated 
with the magnetic field variation rather than the sole direct effect of the magnetic 
field itself. Li, Basu, Sofia, et al. (2003), using models of the structure and evolution 
of the Sun, found that turbulence near the surface of the Sun plays a major role in 
solar variability, and only a model that includes a magnetically modulated turbulent 
mechanism can agree with the observed correlation between the frequency shift and 
the solar cycle. In such a dynamic model, the evolution of the subsurface layers of 
the Sun through the activity cycle plays an important role. 

The frequencies of the global modes give the radius and latitude (r, 0) part of 
the structure. While, local helioseismic techniques such as ring-diagram analysis 
(Hill, 1988) allow the determination of the three-dimensional structure of the Sun, 
allowing the study of localized areas in the solar surface, such as those in active 
regions. Large variations of the mode frequencies observed in and near sunspots 
in comparison to magnetically quiet regions are well known to be correlated with 
variations in the average surface magnetic field between the corresponding regions. 
(e.g., Rajaguru, Basu and Antia, 2001 and Rabello-Soares, Bogart and Basu, 2007). 
Whether the frequencies are changed directly by the magnetic field or indirectly 
through an associated change in the solar structure (such as a pressure change) is 
still a matter of debate. A detailed analysis of the frequency-shift characteristics will 
hopefully help understand their physical origin. 

Basu, Antia and Bogart (2004), using ring-diagram analysis, found that the sound 
speed is lower in the immediate subsurface layers of an active region than of a 
magnetically quiet region, while the opposite is true for depths below about 7 Mm. 
However, Basu and Antia (2002), using global analysis, have found no observable 
structural changes in the inner layers of the Sun below a depth of 21 Mm associated 
with the magnetic-activity induced frequency shifts. They were, however, unable to 
get closer to the solar surface due to the lack of high-degree modes in their mode 
set. High-degree global analysis is important to complete the picture of the near- 
surface layers. Besides, the determination of high-degree frequencies using different 
methods allows us to check the results against each other giving confidence in the 
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results and avoiding systematic errors. We should point out that, although the high- 
degree modes have short lifetimes (one -ten hours for 100 < I < 600 accordingly to 
Burtseva, Kholikov, Serebryanskiy, et at, 2007) propagating only locally, they are 
averaged over most of the solar surface using spherical harmonic decomposition (over 
a relatively long time series) and thus can still be called global analysis. 

In the traditional global helioseismology data-analysis methodology, a time series 
of full-disk Doppler solar images is decomposed into spherical harmonic coefficients, 
characterized by its degree (I) and its azimuthal order (m). Each coefficient time 
series is Fourier transformed, and the order of the radial wave function (n) gets 
separated in the frequency domain. However, a spherical harmonic decomposition is 
not orthonormal over less than the full sphere - i.e., the solar surface that can be 
observed from a single view point- resulting in what is referred as spatial leakage. 
At low and intermediate degrees, most of these leaks are separated in the frequency 
domain from the target mode (except for some m leaks) and individual modes can 
be identified and fitted. However, at high degrees, the spatial leaks lie closer in 
frequency (due to a smaller mode separation) and, at high frequency, the modes 
become wider (as the mode lifetimes get smaller), resulting in the overlap of the 
target mode with the spatial leaks that merges individual peaks into ridges (see 
Figure 1 in Rabello-Soares, Korzennik and Schou, 2001). The characteristics of the 
resulting ridge (central frequency, amplitude, etc. . . ) do not correspond to those of 
the underlying target mode. This has so far hindered the estimation of unbiased 
mode parameters at high degrees. 

To recover the actual mode characteristics, we need a very good estimation of the 
relative amplitude of the spatial leaks present in a given (£, m) power spectrum, also 
known as the leakage matrix, which in turn requires a very good knowledge of the in- 
strumental properties (Rabello-Soares, Korzennik and Schou, 2001). In our previous 
papers (Rabello-Soares, Korzennik and Schou, 2001 and Korzennik, Rabello-Soares 
and Schou, 2004, hereafter KRS), we described in detail the large influence of the 
instrumental properties on the amplitude of the leaks and as a consequence in the 
determination of unbiased high-degree mode parameters. 

In the following, we will first describe the data used in this analysis and the 
ridge-to- mode correction applied to them (Sections 2 and 3). In Section 4, we will 
discuss the influence on the mode parameters of each of the instrumental properties 
that were included in the spherical harmonic decomposition of the solar images. We 
then analyze in Section 5 the characteristics of the high-degree mode frequencies and 
their rotational splittings obtained in this work. Finally, in Section 6, we analyze the 
frequency variation induced by the solar cycle. 

2. Observations and Ridge- Parameter Extraction 

The data used in this work consist of full-disk Dopplergrams obtained at a one- minute 
cadence by the Michelson Doppler Imager (MDI) onboard the Solar and Heliospheric 
Observatory (SOHO). We have used two distinct sets of data. One while MDI was 
operating in its 4" resolution mode, allowing the detection of oscillation modes up 
to I w 1500, which we will call from now on the high-^ data set. This is the so- 
called Dynamics Program observing mode, which is available every year for two to 
three months. The second one (hereafter referred to as the medium-^ data set) using 
MDI Structure Program that provides almost continuous coverage year round. In 



p. 3 



M.C. Rabello-Soares et al. 



Table 1. Details of the analyzed time series. The corresponding relative mean 
values of solar UV spectral irradiance and their standard deviation are also 
listed as an indication of the solar activity level. 



Year Starting Date Duration Solar Index Starting Date 





high-£ set 


high-£ set 


rcl. to max. (in %) 


medium-^ set 


1999 


Mar. 13 


77 days 


40 ± 13 


Feb. 03 


2000 


May 27 


45 days 


69 ±11 


Apr. 10 


2001 


Feb. 28 


90 days 


61 ±14 


Jan. 23 


2002 


Feb. 23 


72 days 


80 ±8 


Mar. 31 


2003 


Oct. 18 


38 days 


51 ±17 


Oct. 28 


2004 


Jul. 04 


65 days 


36 ±11 


Aug. 11 


2005 


Jun. 25 


67 days 


30 ±9 


May 26 



this observing mode, the original Dopplergrams are convolved, onboard the MDI 
instrument, with a Gaussian and subsampled on a 200 x 200 grid, thus reducing the 
telemetry requirements but also limiting the spatial resolution to modes with degree 
£ < 300. 

For the high-£ set, we computed the spherical harmonic decomposition of the MDI 
images for modes with 100 < £ < 900. The resulting time series were Fourier trans- 
formed in small segments (4096 minutes) whose spectra were averaged to produce an 
averaged power spectrumwith a low but adequate frequency resolution to fit the ridge 
while reducing the realization noise. The number of averaged spectra varies with the 
year, from 12 segments in 2003, to 30 in 2001. Most of the known instrumental effects 
relevant to the high-degree analysis were included in the spatial decomposition and 
they will be described in Section 4. The peaks in each {£, to) spectrum were then 
fitted using an asymmetric Lorentzian profile with an additive background term, 
given by ten to the power of a second-degree polynomial in frequency (Equation 5 in 
KRS). Since the number of segments used in the average of each (I, to) spectrum is 
large enough, its y 1 distribution can be approximate by a Gaussian distribution and 
a least-square fitting was used. The fitted Lorentzian profile is characterized by the 
following parameters: frequency (v n t m ): amplitude (A n ^^ m ) , width (^in.i.m) and, 
asymmetry {u n ,i,m)- The asymmetric profile used is equivalent to the one defined 
by Nigam and Kosovichev (1998) where their asymmetric parameter is equal to 
a I (2 — a). The frequency splittings were parametrized in terms of Clebsch-Gordan 
coefficients up to (Ritzwoller and Lavely, 1991). The number of modes analysed 
was reduced without loss of information to the work described in this paper and 
hence was easy to handle. The fitting was carried out only every tenth £ and only for 
some 50 equally spaced to values at each £. The central frequency, i.e. the frequency 
free of splitting effects, is taken to be the frequency given by m = in the splitting 
parametrization. Since the even splitting coefficients are zero on average, it is the 
same as the mean frequency (averaged over to). Notice that the fitted parameters 
are ridge parameters and do not correspond to the associated mode parameters as 
discussed in the introduction. In this study, high-i? time series available from 1999 
until 2005 were used and their properties are listed in Table 1. 

The focus of this work are modes with £ > 300 and therefore it is centered on 
the analysis of the high-i? set. We used however the results obtained by one of us 
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analyzing the medium-/ time series (Schou, 1999) to compare with and complement 
our analysis. Using 72-day long time intervals, the central frequency and splitting 
coefficients for a given (n,£) mode were determined directly by fitting symmetric 
Lorentzian profiles to its power spectra (Schou, 1999). Every p mode up to t « 200 
(up to t « 300 for the / modes) and every m were fitted. The 72-day time series 
that best overlap in time with the high-/ time series were selected and are listed in 
Table 1. Only modes with £ > 20 were used. The mode coverage of both analyses is 
illustrated in Figure 1. 

In the medium-/ analysis, the instrumental effects described in Section 4 were not 
taken into account, and neither were the distortion of the eigenfunctions by the solar 
differential rotation or the horizontal component of the oscillation (both described 
in Section 3). However, in the frequency and degree ranges of the medium-/ analysis, 
the spatial leaks are well separated from the target mode and individual modes can 
be identified and fitted, making the above-mentioned effects not as crucial as they are 
for the high-/ analysis. Another relevant difference between the medium- and high-/ 
analysis is that the medium-/ power spectra were fitted using symmetric profiles, 
which is well known to lead to systematic errors in the frequency measurements 
(e.g. Toutain, Appourchaux, Frohlich, et al, 1998; Basu and Antia, 2000). Larson 
and Schou (2008) have recently re-analyzed the medium-/ time series and reported 
that several of these physical effects result in highly significant changes in the mode 
parameters. Their Figure 1 shows the total correction to be applied to the medium-/ 
frequency and splitting coefficients, a\ and 02, used here. Their Figure 2 illustrates 
the frequency changes due to each of these effects. The first two panels correspond to 
the distortion of the eigenfunctions by the solar differential rotation and the horizon- 
tal component of the leakage matrix respectively. The third, forth, and fifth panels 
correspond to the instrumental effects described here in Sections 4.1, 4.2, and 4.3 
respectively. The panels at the bottom show the difference obtained between fitting 



p. 5 



M.C. Rabello-Soares et al. 



symmetric and asymmetric Lorentzians. Although these corrections are significant, 
they correspond to very small variations in the results presented in this paper and 
do not affect our conclusions, as it will be discussed later (Section 5). 

3. Ridge-to-Mode Correction 

Our methodology to recover the mode characteristics from the ridges observed at 
high-degree and high-frequency power spectra consists in generating and fitting a 
sophisticated model of the underlying modes that contribute to the ridge power 
distribution and deduce the offset (A ) between the ridge properties and the target 
mode (Korzennik, 1998). A synthetic power spectrum is computed for each (£,m) 
mode consisting of several asymmetric Lorentzians for each n: one for the target 
mode (£,m) and one for each spatial leak {£ , m ) with a relative amplitude given 
by the leakage matrix. It has the same frequency resolution as the high-^ set power 
spectra and it is generated for the same set of (£, m) modes. 

We used the complete leakage matrix {i.e., radial and horizontal components) 
where the horizontal-to-vertical displacement ratio is given by Christensen-Dalsgaard 
(2003): GMqL/(RqU}^^), where G is the gravitational constant, M© is the solar 
mass, Rq is the solar radius, lu is the cyclic frequency (u) = 2ttv), and L 2 = £{£+ 1). 
For each {£, m) synthetic power spectrum, we have taken into account the contribu- 
tion of the spurious modes (£ , m) that obey the following expressions: |m — m\ < 10 
and \£' -£\< £ d , where £ d = 12 for £ < 600 and £ d is equal to 0.02 1 rounded to the 
closest integer for £ > 600. Spurious modes with (£' , m) values that differ from those 
of the power spectrum (£, m) by more than the amount specified in the equations 
above have a very small relative amplitude and their contribution to the synthetic 
power spectrum profile is negligible, i.e., their inclusion or not does not affect the 
fitted parameters of the peaks in the power spectrum (Rabello-Soares, Korzennik 
and Schou, 2005). We also included in our model the distortion of the eigenfunctions 
by the solar differential rotation as described in Woodard (1989), using the values 
given by Schou, Antia, Basu, et al. (1998) for the solar rotation. The coefficients in 
this superposition become negligible when \£ — £ \ > £ c where £ c = 10 for £ < 400 
and the next even integer to 3 + 0.02 £ for £ > 400 (Rabello-Soares, Korzennik and 
Schou, 2005). 

The profile resulting from the overlap of several profiles of nearby spatial leaks is 
reasonably well modeled by a single profile when the ratio of the mode width to their 
separation (given by dv/d£) becomes large. The difference is smaller than the rms 
of the observed residuals to the fit, 5% -7% at all degrees, which are dominated by 
the realization noise (KRS). The synthetic power spectrum is then fitted following 
the methodology applied to the high-£ set (Section 2), providing the ridge central 
frequency v^° t dcl and its splittings a™ del (i = 1,6). The offset is given by the differ- 
ence between the modeled ridge and the target mode parameter given by the input 
value of the parameter used to generate the synthetic power spectra. For a given 
mode (n,£), the offsets in central frequency can be then written as: 

A o model input / -i \ 

A v n j = - v n £ . (1) 

The offsets in each of the fitted a-coefficients are obtained in the same manner. 
Realistic input values based on observed mode parameters were used (as described 
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in KRS). The input linewidth is given by the square root of 7m de +7wfj where 7 mo de 
is an estimate of the mode linewidth and 7wf is the width of the window function in 
our ease given by the length of the high-£ set time string (i.e., 4096 minutes). This 
expression was obtained by Korzennik (1990) assuming that the observed power 
spectrum is a convolution of the "true" power spectrum (i.e., the one that would 
have been obtained with a infinite time series), with the power spectrum of the 
window function both represented by Gaussian profiles, which is adequate for the 
purpose at hand. 

If the leakage matrix is correct and complete, our simulations should adequately 
estimate the parameter offsets (A°) and we would be able to obtain corrected mode 
parameters from the observed ridge. Thus in the case of the central frequency, we 
would have: 

corrected observed a /o\ 

v n .i = v n ,i ~ A v n £. (2) 

The corrected mode a-coefficients are obtained as in Equation (2). 

Figure 2 shows the estimated offsets for the central frequency and the splitting 
coefficients 01, ai, and 03. They are primarily a function of frequency. Except for ai, 
the offsets are quite large in comparison with the observed fitting uncertainties of 
the corresponding parameter which, for the central frequency, ai and 613 coefficients, 
are usually smaller than 0.5 /iHz, 1 11Hz, and 0.8 nHz respectively. The width and 
amplitude offsets will be described in a future paper. Accordingly to Korzennik 
(1998), the asymmetry of the ridge seems to be the same as the asymmetry of the 
underlying mode (see also KRS). 

Another methodology to recover the mode characteristics from the observed ridges 
is to fit a sum of Lorentzians or asymmetric profiles to a given ridge: one for the target 
mode and one for each spatial leak that has a significant amplitude (Reiter, Rhodes, 
Kosovichev, et al, 2002; Reiter, Kosovichev, Rhodes, et al, 2003; Reiter, Rhodes, 
Kosovichev, et al., 2004). Such an approach is likely to be particularly appropriate in 
the transition region between well-resolved modes and ridges, where the spatial leaks 
start to overlap with the target mode but do not yet blend fully into a ridge. The 
reasoning behind the method used here is that, when the spatial leaks are completely 
blended into a ridge, there is not enough information in the power spectra to justify 
modeling the individual spatial leaks as part of the fitting of a given observed spec- 
trum. The difference in the profile of a sum of overlapped asymmetric profiles and 
a single asymmetric profile is much smaller than the observed fitting residuals. This 
method is significantly less computationally demanding since the profiles fitted are 
much simpler. As a result we can more easily check its reliability by using different 
leakage matrices or time series produced with a different spatial decomposition, and 
quantitatively validating the corrections; as described in Section 4 and KRS. This 
validation step ensures the reliability of the method and allow us to estimate a 
quantitative upper limit on any residual bias, a crucial step in the intricate high-£ 
analysis. Note that the method described by Reiter et al. also relies on a very good 
estimation of the leakage matrix to obtain unbiased mode parameters. Indeed, the 
same leakage matrices have been used for both methods and it is likely that both 
the random and systematic errors will be quite similar where the modes are fully 
blended into ridges. 
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Figure 2. Estimated offsets for the central frequency (u) and splitting coefficients (oi, «2, and 
C13). The error bars correspond to the uncertainties when fitting the synthetic power spectra. 
Modes with n < 5 are in color with the / modes in red. 



4. Influence of Instrumental Properties on the High-Degree Power 
Spectra 

In the determination of unbiased high-degree mode parameters, it is necessary to have 
a good knowledge of the properties of the instrument used to collect the data. The 
instrumental characteristics must be taken into account either in the image spatial 
decomposition or in the leakage matrix calculation to obtain a correct estimation of 
the amplitude of the spatial leaks. In a continuing effort to infer unbiased estimates 
of high-degree mode parameters using the high-resolution observations from the MDI 
Dynamics Program, the following instrumental properties were included in the spatial 
decomposition of the high-i? data set one at a time and their effect on the observed 
power spectra analyzed: (i) the correct instantaneous image scale, (ii) the radial 
image distortion, (in) the non-radial image distortion, (iv) the effective P angle, and 
(v) a correction to the Carrington elements. Once a given instrumental property is 
analyzed, it is incorporated in the analysis from then on in the paper. For example, 
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in the analysis of the radial image distortion (Section 4.2), the correct image scale 
(Section 4.1) was used and, in the non-radial image distortion analysis (Section 4.3), 
the correct image scale and the radial image distortion were included. 

4.1. Image Scale 

Variations in the amount of defocus of the instrument have a direct influence in the 
image scale at the detector 1 . Although the MDI instrument has been very stable 
over the more than 11 years that it has been observing the Sun, continuous exposure 
to solar radiation has increased the instrument front window absorption resulting in 
a continuous small increase of the instrument defocus. Moreover the change of the 
front-window temperature due to the satellite orbit around the Sun also adds a small 
annual variation in the image defocus. The instrument has however an adjustable 
focus with nine possible positions chosen to best suit various science needs, resulting 
unfortunately in abrupt jumps in the image defocus every time that a new position is 
chosen and which are responsible for the largest variations (see Figure 5 in Rabello- 
Soares, Korzennik and Schou, 2001). The average size change (at the solar limb) 
per focus step is 0.529±0.002 pixels (Kuhn, Bush, Emilio, et al, 2004). Due to 
these different time-varying focus variation, the image scale must be continuously 
estimated and the correct value used in the spatial decomposition. The image scale is 
obtained by measuring the observed image radius, which is defined as the inflection 
point in the radial limb-darkening function. The FNDLMB routine in the GONG 
Reduction and Analysis Software Package (GRASP) was used. It is available from 
the National Solar Observatory, Tucson, AZ, U.S.A.. 

To show the influence of the image scale in the high-£ mode parameters, we 
compared the observed ridge parameters obtained using two different spatial anal- 
yses of the 1999 time series. In one of the analyses, the time-varying image scale 
was obtained for the actual observational period (i.e., the correct image scale) and 
used in the spatial decomposition. In the other one, a constant value obtained from 
observations taken in early 1996 (at the beginning of the mission), i.e., the wrong 
value for the 1999 time series, was used. It is 0.27% larger on average than the actual 
1999 time-varying image scale. A variation in the observed ridge parameter due to 
a change in the data analysis corresponds to an identical variation in its offset (A ) 
in order to obtain the correct mode parameter (see Equation 2). Figure 3 shows the 
corresponding offset variation. The frequency offsets obtained using the correct image 
scale are systematically larger than the ones using a slightly larger image scale and 
their difference increases with frequency (Figure 3). There is no indication of a degree 
dependence. The variation in the frequency offset is larger for the / modes than for p 
modes (respectively upper and lower branches seen at frequencies smaller than 3 mHz 
in Figure 3). This suggests that the image scale correction has a strong effect on the 
horizontal component of the leakage matrix. The a\ coefficients obtained using the 
correct image scale are systematically larger by 1.44 ±0.02 nHz than using an image 
scale 0.27% larger on average. Their difference is plotted in Figure 3 (bottom panel) 
arbitrarily against degree instead of frequency. The effect on the other parameters 
is very small and their mean difference is listed in Table 2. The 2000 data set was 
used to calculate the values in Table 2, except for the variations due to the effect 



Image scale is the ratio between the size of the solar image at the detector and its actual size. 
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of the image scale where the 1999 epoch was used. The values in the table should 
not depend significantly on the observational period used since we are comparing 
variations in the analysis of the same time series. The variations in the frequency 
offset (Figure 3 top) are quite large in comparison with their absolute values shown 
in Figure 2. The theoretical offset estimation in Figure 2 corresponds to the case 
where the correct values for all instrumental effects were taken into account in the 
spherical harmonic decomposition. 

The red points in Figure 3 were calculated taking into account the image scale 
error in the leakage matrix calculation instead of in the image spatial decomposi- 
tion. The leakage matrix is calculated assuming that a constant and 0.27% larger 
image scale than the actual value was used in the image spatial decomposition. The 
variation in the parameter offsets estimated from the synthetic power spectra using 
Equation (1) matchs very well the observations in most cases. For modes above 4 
mHz, the frequency changes were underestimated by this method. This could be 
because we did not take into account the image-scale temporal variation, only the 
average difference for the entire observing period, while the spatial decomposition 
was carried out using a instantaneous image scale estimation. 

Note that the smearing of the image represented by the point-spread function 
(PSF) is not taken into account nor its variation with the amount of defocus. Our 
preliminary analysis indicated that including a simple model of the azimuthally 
averaged estimation of the PSF in the leakage matrix calculation has a very small 
effect on the offsets. It affected mostly the frequency offset and only by less than 
0.2 /iHz (Rabello-Soares, Korzennik and Schou, 2006). Recently we found that the 
observed ridge frequencies obtained for the observational periods where MDI was set 
to a large defocus, i.e., 1996 to 1998, are larger than the ones obtained for the other 
periods where the instrument was nearly in focus, after correcting for the solar-cycle 
variation; their maximum difference is of the order of a few fiRz (Rabello-Soares, 
Korzennik and Schou, 2008). A possible explanation is that an azimuthally averaged 
estimation of the PSF is not a good approximation to the true PSF of the instrument, 
which is known to depend on the azimuth angle (Schou and Bogart, 1998) with a 
phase that changes with focus position (e.g., KRS). Unfortunately, there is not a 
good estimate of the PSF for the MDI Dynamics Program at the moment. We plan 
to use an approximation to the azimuthally-varying PSF and analyze its influence 
in the mode parameter determination in the future. 

4.2. Radial Image Distortion 

The ray-trace model of the MDI optical configuration predicts a radial distortion 
( Ar) which depends on the distance from the CCD center (r) : 

— = b (r 2 - <r> 2 ), (3) 
r 

where b = 7 x I0~ 9 pixels -2 and <r> is the observed image mean radius (Kuhn, 
Bush, Emilio, et al, 2004). The distortion causes the apparent solar radius to be 
larger by wO.17% («0.8 pixels or 17 /jm). Thus, the second term in Equation (3) 
was added to ensure that the distorted and undistorted images have the same mean 
radius. 

As in Figure 3, Figure 4 shows the offset variation corresponding to the difference 
in the observed ridge frequencies when including this distortion in the spatial decom- 
position. Similarly to the image scale, the frequency offset increases with frequency. 
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Figure 3. Difference in the frequency and the ai-coefficient offsets using, in the spatial 
decomposition for the 1999 high-£ set, the correct time- varying image scale ("corr") and a 
constant value, which is 0.27% larger on average than the actual values ("wrong"). The errors 
are given by the fitting uncertainties. The estimation of this effect when it is included in the 
leakage matrix calculation instead (see text for details) is shown in red. 



The similarity is expected, since the radial distortion changes the image scale by 
an amount that is a function of the distance from the CCD center r (Equation 3). 
The difference in the ridge width also increases with frequency from —0.064 ± 0.004 
/iHz at v < 2.5 mHz to 0.28 ± 0.04 fjHz at v > 4.5 mHz, in the sense including 
minus not including the distortion, but it is very small in comparison with the fitting 
uncertainties and barely significant. The effect on the other parameters is small and it 
is listed in Table 2. The ridge modeling of this effect, introducing the radial distortion 
in the leakage matrix, agrees well with the observations (red circles in Figure 4). 

4.3. Non-Radial Image Distortion 

Solar images as observed by MDI have a nearly elliptical shape with a difference 
between the semi-major and the semi-minor axis of about 0.6 pixels. This is consistent 
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Table 2. Mean differences in the observed ridge parameters in the sense 
improved minus unimproved. The weighted average and its standard error 
were calculated over all observed modes (500 in all cases). The mode range 
is shown in Figure 1 (black circles). The weight of each determination was 
taken as inversely proportional to the square of the uncertainty. The un- 
certainties of the central frequency and the splitting coefficients arc given 
by the uncertainties in the Clebsch-Gordan expansion. The uncertainties 
in the width, amplitude, and asymmetry are given by the standard devi- 
ation of the m-average. The first line ("Uncertainty") gives the minimum 
and maximum values of the uncertainties, which are very similar in all 
the analysis presented in the Table. The other lines, from top to bottom, 
are the differences in the image scale, radial image distortion, non-radial 
image distortion, P-anglc and Carrington elements respectively. The am- 
plitude values are given in arbitrary units and the asymmetry parameter 
is dimensionless. 



Parameter 


Frequency 


ai 


a 2 X 10 s 


0-3 




GuHz) 


(nHz) 


(nHz) 


(nHz) 


Uncertainty 


0.05 - 1.4 


0.07 - 2.0 


70 - 20000 


0.2 - 8 


Scale 


4.54 ±0.06" 


1.44 ±0.02 


-25 ±9 


0.03 ± 0.02 


Radial 


0.95 ±0.01° 


0.42 ± 0.02 


-95 ±7 


0.03 ± 0.02 


Non-rad. 


0.133 ± 0.008 


-0.11 ±0.02 


240 ± 10 


0.03 ± 0.01 


P angle 


-0.0012 ±0.006 


0.04 ±0.02 


-0.5 ±8 


0.025 ± 0.016 


Carr. el. 


0.00020 ± 0.005 


-0.018 ±0.02 


-1 ±6 


-0.013 ±0.01 



Parameter 


Width (/iHz) 


Amplitude 


Asymmetry 


Uncertainty 


0.3 - 6 


50 - 10 5 


0.003 - 0.12 


Scale 


0.047 ± 0.007 


-4 ±2 


(-1 ±2) x 10~ 4 


Radial 


-0.018 ± 0.005 a 


0.4 ± 1 


(-7.5 ±0.5) x 10~ 4 


Non-rad. 


-0.044 ± 0.005 


2±1 


(1 ±0.4) x 10~ 4 


P angle 


-0.00026 ± 0.005 


-2 ±2 


(9 ±4) x 10~ 7 


Carr. el. 


-0.0050 ± 0.004 


0.9 ± 1 


(-0.4 ±4) x lO" 6 



a In these cases, the difference depends on the frequency (see text). 



with a small tilt (ss2°) in the detector around an axis that is rotated 56° from the 
detector's horizontal a;-axis and it introduces a non-radial image distortion. In KRS, 
we estimated this distortion using different observational methods. Unfortunately, 
the distortion varies by as much as 35% depending on the data used, with the 
correspondent CCD tilt varying from 1.71° to 2.6°. Kuhn, Bush, Emilio, et at (2004) 
also estimated the non-radial distortion using yet another observational method - 
the Mercury transit across the Sun on 7 May 2003 — and found a 3.3° tilt. Although 
their equations to estimate the distortion have the same general form as ours (KRS), 
there are other differences besides the tilt angle between the two calculations. Their 
estimation of the distortion varies by 14% in relation to our estimation using a 2.6° 
tilt. These discrepancies in the estimated distortion using different observational 
methods might be attributed to a number of reasons such as the inaccuracy of our 
simple model for the non-circular shape of the solar images or the influence of an 
optical aberration (an asymmetric PSF, for example). 
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To analyze the effect of the non-radial distortion on the power spectra, we included 
in the spherical harmonic decomposition our estimation that has the largest tilt angle 
(2.6°), which corresponds to the distortion that better reproduces the solar limb 
shape. Although, to the moment, we were unable to determine precisely the non- 
radial distortion pattern, our estimation provides an improvement to the analysis 
and it will be incorporated from then on in the paper. Fortunately, it has an overall 
small effect on the observed spectra, as can be inferred from Table 2, and we can safely 
extrapolate that small variations from this distortion pattern can only correspond to 
variations in the parameters smaller than the differences shown in the table which 
were obtained comparing with using no correction for the non-radial distortion, and 
most likely negligible. 

4.4. Position Angle P 

The roll angle of the SOHO spacecraft is maintained such that the effective position 
angle 2 (P off ) of the MDI images should always be zero. However, a 0.2° difference 
has been measured by intereomparing MDI and GONG images obtained in 1999 
and 2000 (Cliff Toner, private communication). Beck and Giles (2005) estimated a 
0.07° difference, assuming that there is no equator-crossing flow, after re-gaining 
contact with the SOHO spacecraft (in 1999) and 0.1° before losting contact. No 
noticeable effect is seen in the observed ridge parameters after including a 0.2° 
correction (Table 2). We do not see the 3 nHz variation in the ai-coefficient offset that 
was predicted by our ridge model using a slightly higher correction of 0.25° (KRS). 
This is probably because we did not accurately model the _P-angle correction in the 
leakage matrix calculation, but added a very simple approximation of its effect. 



2 It is the position angle of the northern extremity of the solar rotation axis with respect to 
MDI detector y-axis. 
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4.5. Carrington Elements 

Accordingly to Beck and Giles (2005), the standard values used for the two angles 
specifying the orientation of the solar rotation axis (i, f2) 3 , known as the Carrington 
elements, are off by Ai = 0.095° ± 0.002° and Afl = 0.17° ±0.1°. This introduces a 
time- varying correction in the calculation of the rotation axis projection in the plane 
of the sky, i.e. the position angle, P, and the roll angle, -Bo 4 . This correction in the 
P angle will be on top of the one mentioned in the previous section. Introducing a 
correction of Ai = 0.1° and Af2 = 0.1° in the image spatial decomposition has no 
significant effect on the observed ridge parameters (Table 2). 

5. The High-Degree Mode Parameters 

Here we analyze the frequencies and splitting coefficients obtained using the high-£ 
data set. The mode parameters were obtained after including the five known in- 
strumental effects in the spatial decomposition of the high-£ data sets described in 
Section 4, fitting their observed power spectra using an asymmetric Lorentzian profile 
(Section 2) and applying the ridge-to-mode correction to the fitted ridge parameters 
(Section 3). 

First, in Section 5.1, the high-£ data set mode frequencies and splitting coefficients 
are compared with the values obtained by the medium-^ analysis. Figure 1 shows the 
region of common modes between the two sets. This comparison is done to check the 
goodness of the estimation of the mode parameters from the observed ridge in the 
high-£ data set at these high medium-^ common modes. Then, the estimated high-£ 
frequencies and splitting coefficients at £ > 100 are analyzed (Sections 5.2 and 5.3). 

5.1. Medium- and High-£ Set Comparison 

Figure 5 shows the differences between the mode frequency and the splitting co- 
efficients obtained using the 72-day medium-^ data set described in Section 2 and 
using the high-£ set both observed during 2004 (Table 1). The high-£ set used in this 
comparison was analyzed as described in Sections 2, 3, and including all instrumental 
effects analyzed in the Section 4. In the medium-^ mode range, the spatial leaks are 
well separated from the target mode and individual modes can be identified and 
fitted. By decreasing the observed high-£ set frequency resolution using a short time 
string (4096 minutes), and thus increasing the width of the window function, we 
force the width of the spatial leaks to increase. The spatial leaks in the high-£ power 
spectra now overlap with the target mode forming a ridge at £ as low as 100 and the 
ridge-to-mode correction described in Section 3 is applied. To increase the number 
of common modes, the high-€ power spectra used in the comparison were fitted for 
every £ (and not every tenth £). The set of modes used in this comparison consists 
of 420 p modes with 100 < £ < 200 and 1.7 < v < 3.5 mHz and 60 /-modes with 
230 < I < 300 and 1.5 < v < 1.8 mHz observed during 2004. 



3 j is the angle between the plane of the ecliptic and the solar Equator and J7 is the angle 
between the crossing point of the solar Equator with the ecliptic and the Vernal Equinox. 

4 Bo is the heliographic latitude of the central point of the disk and presents an annual variation. 

5 An even shorter time series, 2048 minutes, was used to check the results at the lowest high-£ 
modes. 
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Figure 5. Central frequency and splitting coefficients difference between the medium- and 
high-4? sets observed during 2004 (Table 1). The differences are in the sense high- minus 
medium-^ set. Modes with n < 5 are in color where the / modes are in red. 

The mode frequencies agree within 1 /iHz (Figure 5 top panel). Their difference 
varies with frequency having a maximum at 2 mHz and decreasing to zero (or near 
zero) at 1.6 and 3 mHz. There is no indication of a degree dependence. Although 
1 /iHz is a small value, it is almost one order of magnitude larger than the fitting 
uncertainty of these modes (0.12 /iHz on average) and hence undesirable. This seems 
to indicate that the frequency offsets /S.°v (Figure 2) are underestimated by «35% 
for modes with frequencies in the range 1.7 to 2.3 mHz, this amount decreases to 
10% around 1.6 mHz and zero in the interval between 2.6 and 3.2 mHz where the 
estimation is correct. There is some indication that at 3.4 mHz the offsets are again 
underestimated (by ~30%). Unfortunately, there are no common modes at higher 
frequencies. 

A variation in the mode frequency that is purely a function of frequency does not 
affect the outcome of the solar structure inversion, since it is removed together with 
the near-surface errors in the physics of the solar model (see Section 5.2). The instru- 
mental effect that is probably causing the frequency offset to be underestimated at 
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certain frequencies might be the instrumental PSF which is not included in the image 
spatial decomposition or in the leakage matrix calculation (see Section 4.1). Besides 
the frequency dependence, there is a puzzling frequency difference (« 0.5/iHz) in the 
P2 modes (green circles in Figure 5 top panel) with respect to the adjacent n values. 

The difference in the a\ coefficients between the medium and high-£ set is negligi- 
ble except for n = 0, 1 modes (second panel in Figure 5). Their mean difference 
normalized by the fitting uncertainties is, in units of er, —0.2 for n > 2, 2 for 
n = 1, and 3 for the / modes. The horizontal-to- vertical displacement ratio decreases 
exponentially with n. Thus, the differences observed only for n = 0, 1 modes could 
be an indication of an unaccounted instrumental effect that has a stronger effect 
on the horizontal component of the leakage matrix than on its vertical component. 
Note also that the / modes in the medium-^ set are fitted using a different frequency 
interval around the peak in the power spectrum than the p modes (Schou, 1999). 

There are small differences between the medium- and high-^ set for the splitting 
coefficients a 2 and 03 (third and bottom panels in Figure 5). Their mean differences 
normalized by the fitting uncertainties of the high-^ power spectra are ai = 1 and 
03 = —1 in units of a. 

In conclusion, the estimation of the mode frequency and splitting coefficients 
from the observed ridge in the high-^ data set at these moderate-degree values is, in 
general, quite good. However, there is still room for improvement, specially for the 
central frequency and the /-mode a% splitting coefficient. 

The corrections in the medium-^ mode parameters calculated by Larson and Schou 
(2008), described in Section 2, do not have a large influence in our results. Their 
effect in Figure 5, including the variation between fitting symmetric and asymmetric 
profiles, is small in comparison with the difference between the medium- and high-/! 
sets. The mean total correction in the medium-^ frequency and splitting coefficients 
ai, a2 and 03 to be applied to the values used here are: 0.130 ±0.004 /iHz, - 
0.244 ± 0.009 nHz, (1.78 ± 0.08) x 10~ 3 nHz and -0.137 ± 0.006 nHz respectively (Lar- 
son and Schou, 2008). The average was calculated over the 420 medium- and high-^ 
common modes. The largest correction is in the central frequency. The mean fre- 
quency variation between fitting symmetric and asymmetric profiles is 0.043 ± 0.002 
/iHz for the common modes (Larson and Schou, 2008). The improved medium-^ 
analysis will decrease slightly the differences between medium- and high-^ frequencies 
showed in the top panel of Figure 5. However, it does not change the overall behavior 
of the frequency differences. 

5.2. The Central Frequency 

Figure 6 shows the difference between the high-^ set mode frequencies (100 < I < 
900) and their theoretical value calculated from Christensen-Dalsgaard's model S 
(Christensen-Dalsgaard, Dappen, Ajukov, et at, 1996) as a function of degree for 
different n values. The high-£ set used in this comparison was analyzed as described 
in Sections 2, 3, and including all instrumental effects analyzed in the Section 4. The 
mode range is shown in Figure 1. The medium-^ set frequencies are also plotted as 
a reference. In the absence of any acceptable theory to describe the physics of the 
layers near the solar photosphere, the difference between the observed and theoretical 
frequencies due to the near-surface errors in the model are well known to be quite 
large {e.g., Christensen-Dalsgaard and Berthomieu, 1991). The general trend is such 
that the observed p-mode frequencies are smaller than their theoretical prediction 
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and, at a high enough degree, the differences increase with degree and with frequency. 
Accordingly to the results obtained by the high-£ set, this difference can be as large 
as 60 /iHz. In fact, for a given n, the frequency differences increase almost linearly 
with degree with a slope that also increases linearly with n. The high-^ set analysis 
gives consistent frequencies for all observational periods listed in Table 1, where the 
differences in the frequencies obtained at the different observational periods can be 
explained by their well known variation with solar cycle activity and it is described 
in detail in Section 6. 

In order to compensate for the frequency shifts due to the near-surface errors 
in the model, an unknown function (-F sur f), the so-called surface term (F aul {) is 
usually added to the equation governing helioseismic inversions: F SUI f = Qn.e x 
{^ h 'n,l/ l/ n,l)suif, where Q n ,e is the mode inertia normalized by the inertia of a radial 
mode of the same frequency. Christensen-Dalsgaard, Thompson and Gough (1989) 
from the asymptotic theory of solar p modes pointed out that low- and moderate- 
degree modes propagate nearly vertically near the surface; thus, their behavior in 
this region is essentially independent of degree and depends only on frequency. This 
however does not hold for high-degree modes as can be inferred from Figure 6. A 
second-order asymptotic approximation, where the surface term is a function not 
only of frequency but also of L/v n g (Brodsky and Vorontsov, 1993), must be used 
as shown by Di Mauro, Christensen-Dalsgaard, Rabello-Soares, et al. (2002). 

Figure 6 also shows for comparison the frequency difference applying the ring- 
analysis technique obtained by Rabello-Soares, Bogart and Basu (2007). A magnet- 
ically quiet region (15° in diameter) observed by the GONG+ network on 4 July 
2005 during Carrington rotation 2031 centered at 115o in longitude and 3o south 
in latitude was tracked for 8192 minutes at the appropriate photospheric rotation 
rate. The region crossed the central meridian at the middle of its tracking interval. 
The tracked region was mapped to a plane using Postel's projection and its power 
spectra, given by the three-dimensional Fourier transform of the temporal series of 
images, were fitted using the 13-parameter model of Basu and Antia (1999). The 
wavenumber (k) can be identified with the degree of a spherical harmonic mode 
of global oscillations by L = fcJ?©. As the oscillations in a plane-parallel geometry 
are only discrete in radial order, I does not need to be an integer. Each "mode" 
is obtained by fitting a region of power spectrum that has significant overlap with 
those covered by neighboring "modes", making them not strictly independent. To 
confirm that the frequencies obtained for this particular solar region correspond to 
typical values, they were compared with the frequencies obtained for 13 additional 
quiet regions, with latitudes ranging from -12° to +12° (Rabello-Soares, Bogart 
and Basu, 2007). Their difference is smaller than 3 ^iHz for n < 5. Part of this 
variation is due to the fact that the projection of the spherical solar surface onto 
a flat detector introduces some foreshortening that depends on the distance of the 
region from disk center, which can introduce systematic errors in determining the 
mode characteristics. 

The frequency differences (in relation to the theoretical value) obtained using ring 
analysis present the same trend as using spherical harmonic decomposition, except 
for the / modes. In most cases, the ring-analysis frequencies are smaller than the 
ones obtained by the high-i? set. For p modes, their difference in the sense global 
minus ring analysis varies between -4 and 6 ^iHz for n < 5 and it could be as large 
as 20 /iHz for n > 5. Due to the small size of the region analyzed, the ring-analysis 
power spectra have low spatial resolution doing poorly at medium degree and at 
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Figure 6. Difference between observed frequencies and their theoretical value as a function 
of degree for different n values using the high-£ set data obtained in 2004 (in black). The 
differences are in the sense observed minus theoretical. For comparison, the results for the 
medium-^ set (red crosses) and from ring-diagram analysis (green stars) are also plotted. The 
fitting uncertainties are given by the error bars. Note the different scale for the / modes. 



high n, as can be seen for n = 6, 7 in Figure 6. Note also that for n = 4 and 
I < 300 modes in Figure 6, the difference between the ring-analysis and theoretical 
frequencies is constant and does not continue to decrease as £ decreases. This is 
probably an indication of a poor frequency determination by the ring analysis at 
these medium-ii values. In conclusion, for the p modes, the ring-analysis frequencies 
agree with the global-analysis values within ±6 /iHz, which is much smaller than the 
difference between either of them and the frequencies obtained from the solar model. 

For the / modes, the frequency differences (in relation to the theoretical value) ob- 
tained using ring and global analysis are different. The medium-^ /-mode frequency 
is on average 0.9 ±0.1 /iHz larger than the theoretical values. The /-mode frequency 
obtained by ring analysis for £ < 450 is on average 1.1 ± 0.6 /iHz larger than the 
theoretical frequency. For I > 450, it decreases sharply with degree and it can be 
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described by a fifth degree polynomial . At £ = 880, it is 13 pHz smaller than the 
theoretical values. Such a sharp decrease in the /-mode frequency in relation to the 
theoretical values similar to the one observed for the p-modes was already reported 
by other authors (Duvall, Kosovichev and Murawski, 1998 and references within). 
The frequencies of the / modes, contrary to the p modes, depend only weakly on 
the hydrostatic structure of the model (e.g. Gough, 1993) and a possible explanation 
is that the /-mode frequencies are reduced by granulation (Murawski and Roberts, 
1993). However, the /-mode frequency obtained by global analysis presents a different 
trend in relation to the theoretical frequency. For £ < 71 0, it is only 2.7 ± 0.5 /zHz on 
average larger than the theoretical frequency. At larger £, it decreases linearly with 
degree with a slope of —0.021 ±0.001 /iHz: it is zero at £ = 850. The global-analysis 
results suggest a much weaker interaction between the / modes and the granulation 
than predicted by the ring analysis. 

The error bars in Figure 6, given by the fitting uncertainties are too small to 
be seen, except for the ring-analysis frequencies at n > 6. The observed frequency 
errors obtained using the medium-^ and high-^! sets and ring analysis are in the range 
0.007-0.4 /iHz, 0.06-0.6 /iHz and 0.2-15 fiBz respectively. The ratio between the 
global- and ring-analysis fitting uncertainties varies from 0.9 to 33. The ring-analysis 
uncertainties are similar to those of global analysis at low n (n < 2) and high £ 
(£ > 700). They become much higher than the global-analysis uncertainties as n 
increases and as £ decreases. The frequency is not a parameter in the model used 
to fit the ring-analysis power spectra (see Basu and Antia, 1999), the observational 
error estimated for the fitted width is used instead, following Rajaguru, Basu and 
Antia (2001). 

5.3. The Rotational Splitting Coefficients 

High-i splittings can be used to infer the solar rotation rate in the outermost lay- 
ers of the solar convection zone. At these layers, there is a steep gradient of the 
rotation, making it a very interesting region to study. This was first suggested by 
the fact that different surface markers such as sunspots, faculae, Ha filaments, and 
supergranular network present different rotation rates [e.g. Snodgrass, 1992) which 
has been interpreted by assuming that the markers are anchored at different depths 
(Foukal, 1972). The solar rotation rate determined by helioseismology (using modes 
with £ < 300) has a local maximum at about O.95i?0 (35 Mm below the solar 
surface) and decreases fast towards the surface (Howe, Christensen-Dalsgaard, Hill, 
et al., 2007). Using /-mode splittings with degree between 117 and 300, Corbard and 
Thompson (2002) estimated a constant radial gradient of solar angular velocity of 
around —400 nHz/i?0 in the outer 15 Mm between the Equator and 30° latitude. 
Above 30° it decreases in absolute magnitude to zero around 50°. The addition of 
high-degree splittings to this analysis will help to constrain the determination of the 
solar rotation profile both closer to the solar surface and at higher latitudes. 

Figure 7 shows a smooth variation with degree of the high-i! set a\ coefficients 
as expected, implying a good estimation of the mode splittings. For a given n, 
the p-mode a\ coefficients decrease with degree until £ rj 300 by 1 - 2% (which 



6 The polynomial in I has the following coefficients in units of /iHz: cq = 350 ± 80, c\ = 
-2.7 ± 0.7, c 2 = 0.008 ± 0.002, c 3 = (-1.3 ± 0.3) x 10~ 5 , c 4 = (1 ± 0.2) x 10" 8 and c 5 = 
(-3.2 ± 0.7) x 10" 12 . 
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corresponds to 4-6 nHz) where it becomes approximately constant. And the /- 
mode ai-coefficients decrease almost linearly with degree from £ w 300 to 900 by 
2.5% (or 9 nHz). For a given n, the a\ variation with £ shown in Figure 7 presents the 
same trend if plotted against the location of the mode lower turning point (given by 
L/u), indicating the expected decrease of the solar rotation with solar radius near 
the solar surface. The p-mode ai-coefficients, especially for 2 < n < 5, present a 
sharp decrease near £ = 200 correspoding to a turning point of depth of «30 Mm, 
which is at the beginning of the solar rotation decrease. 

Although the medium-^ a3-coefncient decreases slightly in absolute value with 
degree (from -7.8 nHz at £ w 50 to -7.3 nHz at £ « 300), the high-l? a3-coefficient 
mean value (-8.6 ± 1.2 nHz) is larger in absolute value than the mean medium-^ value 
(-7.5 ±0.3 nHz). The high-^ a3-coefficient does not show any variation with degree. 
However, it has a small absolute value at low frequency [y < 2.2 mHz), which is 
similar to the medium-^ mean: -7.4±0.5 nHz. 

While the odd coefficients provide information about the internal solar rotation 
rate, the even coefficients arise from latitudinal structural variation, centrifugal dis- 
tortion, and magnetic fields. The (^-coefficients obtained with the high-£ set are on 
average zero (1500 < v < 5200 /iHz). The a2-coefficients obtained with the medium- 
£ set are on average zero for v < 2 mHz. At higher frequencies, they start to decrease 
with frequency, where for v = 3550 ± 50 /iHz, its mean value is -0.08 ±0.03 nHz. 

6. Frequency Variation Over the Solar Cycle 

The correlation between solar acoustic mode frequencies and the magnetic activity 
cycle is well established and has been substantially studied during the last and cur- 
rent solar cycles (see for example Chaplin, Elsworth, Miller, et al, 2007, Dziembowski 
and Goode, 2005 and references therein). However, its physical origin is still a matter 
of debate and the detailed analysis of the frequency shift characteristics are likely 
to contribute to solving this problem. The data sets used in this work and listed 
in Table 1 cover a considerable part of solar cycle 23. They also cover a large 
degree range (20 < £ < 900) when combining the medium- and high-£ sets. In a 
previous work (Rabello-Soares, Korzennik and Schou, 2008) , we used these two data 
sets to analyze the characteristics of the frequency variation along the solar cycle. 
Here we improved our previous analysis by calculating the probability that a linear 
relationship indeed exists between the frequency shift of a given (n, £) mode and the 
solar-activity index. Although the medium- and high-£ sets are fitted with symmetric 
and asymmetric profiles respectively, the resulting difference in frequency does not 
change on average over the solar cycle (Larson and Schou, 2008). 

The frequency variation with the solar cycle is given by the frequency difference 
between 1999-2004 epochs with respect to the 2005 epoch, the one with the lowest 
activity index in our data set. For each (n, £) mode, their frequency shifts were fitted 
assuming a linear relationship with the solar activity index with a zero intercept 
and using a weighted least-squares minimization. For a detailed description, see 
Rabello-Soares, Korzennik and Schou (2008). Four different solar activity indices 
commonly used in the literature and that are available in the period 1999-2005 
were used in this analysis: the solar UV spectral irradiance (given by the NOAA 



7 http: / /www. ngdc.noaa.gov/stp/SOLAR/ftpsolaruv. htm I 
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Figure 7. The a\ rotational-splitting coefficients as a function of degree for the 2004 data 
set, where the ones obtained with the medium-£ set are in red and with the high-£ set in black. 
The fitting uncertainties are given by the error bars. 



Mg II core-to-wing ratio), the Magnetic Plage Strength Index (MPSI 8 , Mt. Wilson 
Observatory), the solar-radio 10.7-cm flux 9 (National Research Council of Canada) 
and the sunspot number (SSN 10 , SIDC, RWC, Bel gium) . The minimum-to-maximum 
solar cycle frequency shift (6i/^ e ) was estimated by multiplying the slope of the linear 
fit by the corresponding solar index variation between the maximum and minimum of 
cycle 23. The mean activity level for each epoch in relation to the maximum reached 
in January 2002 is listed in Table 1, given by the solar UV spectral irradiance. 

In an attempt to include in the analysis of the solar-cycle induced frequency 
shifts only modes that are, in fact, correlated with the solar activity, Rabello-Soares, 
Korzennik and Schou (2008) rejected modes whose Pearson correlation coefficient 



8 http: / /www. astro. ucla.edu/~obs/intro. html 

9 http: / /www. ngdc.noaa.gov/stp/SOLAR/ftpsolarradio. html 

1() http:/ /side, oma.be 
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Figure 8. Probability P u (n,£;i) as a function of the mode frequency for the medium (red) 
and high-^ 1 (black) sets, using the sunspot number as the solar activity index i. 

is smaller than 0.8, or whose slope uncertainty is larger than 20% of its absolute 
value. However, the correlation coefficient cannot be used directly to indicate the 
probability that a linear relationship exists between two observed quantities, i.e., 
if, indeed, there exists a physical linear relation between them. A small correlation 
coefficient might indicate only a small slope in their linear relation. A better approach 
is to calculate the probability that the data points represent a sample derived from an 
uncorrelated parent population. We estimated the probability P u (f) that a random 
sample of N uncorrelated data points would yield a linear-correlation coefficient as 
large as or larger than the observed absolute value of r (Bevington and Robinson, 
1992). If P u (f) is very small, it would indicate that it is very improbable that they 
are linearly uncorrelated. Thus, the probability is high that the frequency shift of a 
particular (n, £) mode and the solar-activity index used are correlated and the linear 
fit is justified. 

For each (n, £) mode and solar-activity index i, we calculated the Pearson cor- 
relation coefficient r(n,£;i) and the probability P u {r) = P u (n,£;i). Figure 8 shows 
the probability P u (n,£;i) as function of the mode frequency. It is very similar for 
all four solar-activity indices. Modes with frequency around 3 mHz have the highest 
probability of being correlated. It decreases for modes with frequency smaller than 
f»2.5 mHz or larger than f»4.5 mHz. There is no noticeable dependence of P u with 
£ or v/£. 

Libbrecht and Woodard (1990) were the first ones to suggest that the solar- 
cycle frequency variation (5v n j) is linearly proportional to the inverse mode inertia 
(I n ,e) using observations obtained at the Big Bear Solar Observatory. The observed 
frequency shift is larger for higher-frequency and for higher-^ modes. As pointed 
out by Libbrecht and Woodard (1990), these modes are more sensitive to surface 
perturbations, because they have, respectively, higher upper- and lower-reflection 
points in the Sun. This indicates that the dominant structural changes during the 
solar cycle, inasmuch as they affect p-mode frequencies, occur near the solar sur- 
face. Chaplin, Appourchaux, Elsworth, et al. (2001) obtained a similar shift at low 
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frequencies (v < 2.5 mHz): 5v n _ ( oc 1/I„j. However, at higher frequencies, they 
found: &v n g oc v a /I n j, where a = 1.91 ± 0.03. They used data sets with a much 
higher duty cycle obtained by the ground-based networks, GONG and BiSON, and 
the SOHO satellite (VIRGO/LOI). 

The value of a depends upon the physical mechanism responsible for affecting 
the acoustic mode frequencies during the solar cycle. Several authors (Gough, 1990; 
Libbrecht and Woodard, 1990; Goldreich, Murray, Willette, et at, 1991) estimated 
a = 3 assuming modifications in the thermal structure of the Sun in a thin layer at 
the photosphere. On the other hand, Gough (1990) pointed out that if the frequency 
shifts are caused by variations in the efficacy of the convection during the solar cycle, 
a « -1. 

Rabello-Soares, Korzennik and Schou (2008) showed that the scaled frequency 
shift can be described with a simple power law at all frequencies: 

re n { v n,lV ,.s 
OV n ,i = Cj —= , (4) 

where Q n ^ is the mode inertia, normalized by the inertia of a radial mode of the 
same frequency, I(,=o(y n g), calculated from Christensen-Dalsgaard's model S (see 
Christensen-Dalsgaard, Dappen, Ajukov, et at, 1996). The / andp modes were fitted 
independently of one another. Chaplin, Appourchaux, Elsworth, et al. (2001) also 
plotted the frequency shift scaled by the normalized mode inertia (upper right-hand 
panel of their Figure 1). However, they chose to fit the frequency shift scaled by the 
mode inertia instead. The multiplication of the frequency shift by Q n ^ normalizes 
the shift to its expected radial equivalent. 

As in Rabello-Soares, Korzennik and Schou (2008), Equation (4) was fitted to 
estimate 7 P (for the p modes) and 7/ (/ modes), using a weighted least-squares 
minimization (thin green line and black dashed line in Figure 9 respectively). Only 
modes that have a probability of 0.1% or less that its variation is linearly uncorrelated 
with the solar index were included in the present analysis. They represent 80% of the 
high-/ set modes and 65% of the medium-/ set modes. The frequency ranges used 
in the fitting of the medium- and high-/ p modes and high-/ / modes are: 2.5-4.1 
mHz, 2.8-4.9 mHz and 1.8-3.0 mHz respectively. It seems that there is a step in 
the p-mode frequency shift around 2.3 mHz and only modes with v > 2.5 mHz were 
included in the fitting of the medium-/ set. p-modes with v < 2.3 mHz seem to 
have a similar slope to the fitted high-frequency modes, but a different y-intercept. 
This is illustrated by the thick green short line in Figure 9. The fact that the / 
modes are affected by the solar cycle in a different way than the p modes is expected 
since they have very different properties. The /-mode is essentially a surface wave 
and its frequency is less likely to be influenced by the solar stratification than the 
p-mode frequency. We repeated the analysis for all four solar indices. The results 
are plotted in Figure 10. The exponents obtained using the different solar indices 
agree within their fitting uncertainty. The weighted mean from the four solar indices 
is 7/ = 1.29 ± 0.07 for the / modes and 7p = 3.60 ± 0.01 for the p modes. 7 P 
was calculated after averaging the values obtained for the medium- and high-/ sets 
(given by the circles in Figure 10). Only / modes obtained with the high-/ set were 
used to estimate 7^. There are only a few modes (seven) in the medium-/ set with 
a probability of 0.1% or less of being uncorrelated, they have a very small frequency 
shift (0.05 ±0.01 fiHz on average), and by consequence they were excluded from the 
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Figure 9. Frequency shift multiplied by the normalized mode inertia as a function of frequency 
obtained, in this case, using the sunspot number as a proxy of solar activity. The long green 
continuous line is the fit to the all p modes with v > 2.5 mHz and the dashed black line is to 
the / modes in the high-£ set. The thick green short line is the fit to p modes with v < 2.3 
mHz using the same slope as the long green line but a different {/-intercept. The blue line 
corresponds to a frequency shift given by a = with an arbitrarily chosen y-intercept. 

fit. The / modes in the medium-^ set have very small frequencies (y < 1.5 mHz) 
and, accordingly to Figure 8, low-frequency modes are not well correlated with the 
solar cycle. High-^ /-modes with v < 1.8 mHz are also not well correlated with the 
solar cycle and were excluded from the analysis. Their frequency shifts show a steep 
increase from values similar to those of the well-correlated medium-£ /-modes at 1.7 
mHz to values similar to the diamonds in Figure 9 at 1.8 mHz. 

In Rabello-Soares, Korzennik and Schou (2008) only the solar UV spectral ir- 
radiance was used as a proxy for the solar-cycle index to calculate the exponent 
7. The value obtained here for the p modes, using the same solar index but a 
better criteria for selecting the modes included in the analysis, is the same as before: 
7 P = 3.63 ± 0.02. However, the value obtained for the / modes (7/ = 1.33 ± 0.2) is 
20% smaller than our previous result. The distinct mode selection criteria account 
for this difference. As the number of / modes is small (45), the fitting of jf is more 
sensitive to the mode selection. 

To compare the above mentioned values of a obtained by Chaplin, Appourchaux, 
Elsworth, et al. (2001) with our results, v" g was divided by Ii = o(v ni i) and then 
fitted by Equation (4) to estimate the corresponding 7 exponent. The a values 
obtained at the two frequency intervals 1.6 < v < 2.5 mHz [a = 0) and 2.5 < 
v < 3.9 mHz (a = 1.91) correspond to 7 = 7.59 ± 0.18 and 7 = 3.58 ± 0.03 
respectively. Accordingly to the authors, the division of the data into two parts 
has been made on a purely "artificial" basis. Rabello-Soares, Korzennik and Schou 
(2008) showed that this imposed "breakpoint" can be explained by the fact that 
\og{5v^ L i I n ,t) has a quadratic dependence on log(i/ n ,f ), with an inflection point at 
2.59 mHz (see Figure 6 in Rabello-Soares, Korzennik and Schou, 2008). Chaplin, 
Appourchaux, Elsworth, et al. (2001) used 10.7-cm radio flux as the solar activity 
index. Their results at the high frequency interval agrees well with ours (represented 
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Figure 10. Exponent 7 obtained fitting Equation (4) to the / (top) and p (bottom) modes 
for four different solar indices. In the top panel, the exponents were obtained using only the 
high-£ set. In the bottom panel, the diamonds and stars were obtained using the medium- and 
high-£ sets respectively and their weighted average is given by the circles. The square is the 7 
corresponding to the a value quoted in Chaplin et al. (2001) for the high-frequency interval. 
The fitting uncertainties arc given by the error bars. The dashed lines are the weighted average 
between the four solar indices. 

by a square in Figure 10 bottom panel). However, at low frequency, their result 
is twice as large as the one obtained here (blue line in Figure 9). The difference 
could be due to the mode selection used to obtain our results. Dziembowski and 
Goode (2005) used also a simple frequency difference including all modes. They 
fitted Sv^ e I n j using truncated Legendre polynomial series and fitted the / and p 
modes independently of one another. For p modes, their fitting agrees with Chaplin, 
Appourchaux, Elsworth, et al. (2001) thus disagreeing with ours at low frequency. 
For / modes, in the frequency range 1.37- 1.60 mHz, the corresponding 7/ (inferring 
from their Fig. 2) is five times large as our determination at v > 1.8 mHz, and it is 
similar to the low-frequency 7 P obtained by Chaplin, Appourchaux, Elsworth, et al. 
(2001). Dziembowski and Goode (2005) noted that the frequency shifts normalized by 
I n l present an opposite trend for the / and p modes. They increase with increasing 
frequency for p-modes and decrease for / modes. This opposite trend disappears 
when the frequency shifts are normalized by Q ni i instead. 

Several authors have reported a sharp decrease of the frequency shifts at high 
frequency. The positive shifts suddenly drop to zero, and become negative reaching 
absolute values much larger than the positive shifts at moderate frequency. The falloff 
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was reported to happen around 3.7 mHz by Jefferies (1998) - using 100 < £ < 250 
modes - and Salabert, Fossat, Gelly, et al. (2004) - using £ < 3 modes, which is 
supported by Libbrecht and Woodard (1990) determinations (5 < £ < 60). However, 
Howe et al. (2002) do not see any falloff for v < 4 mHz modes (£ < 300) and Rhodes, 
Reiter and Schou (2002) see the drop at 5 mHz (using £ < 1000 modes). The exact 
frequencies where the frequency shifts become zero or negative or reach a maximum 
negative value also varies widely between the publications. This sharp decrease in the 
frequency shift has been associated to an increase in the chromospheric temperature 
(e.g., Goldreich, Murray, Willette, et al, 1991 and Jain and Roberts, 1996). From 
Figure 8, high-frequency modes have a large probability P u of being uncorrelated 
with the solar cycle, at least until 5.2 mHz. Note that Pu(r) is calculated for the 
absolute value of r, estimating the probability of being correlated or anti-correlated 
with solar cycle. In our analysis, the frequency shift drops sharply around 4.6 mHz 
and there are a few modes (fs ten) with a negative frequency shift with frequencies 
between 5 and 5.35 mHz. However, only modes with frequency smaller than 4.9 mHz 
obey the mode selection criteria (i.e., P u < 0.1%) and were included in Figure 9. 
The frequency shift drop from 4.6 to 4.9 mHz is not seen in Figure 9 because of the 
mode-inertia normalization. 

As first suggested by Libbrecht and Woodard (1990), the most significant sources 
of frequency shift must be localized near the solar surface since the frequency shifts 
seem to be independent of £. Consider a small perturbation in the solar equilibrium 
model localized near the solar surface and the corresponding changes in the mode 
frequency, 5v n ^. For modes extending substantially more deeply than the region of 
the perturbation, the eigenfunctions are nearly independent of £ at fixed frequency 
in that region (see Figure 8 and the associated discussion in Christensen-Dalsgaard 
and Berthomieu, 1991). From this, it can be inferred that if the quantity 5v n ^Qn,t is 
independent of £ at fixed v for a given set of modes, then the perturbation is probably 
largely localized outside the radius given by the maximum lower turning point of the 
set of modes considered (Christensen-Dalsgaard and Berthomieu, 1991). In the case 
of the frequency shift induced by the solar cycle, Si/ n gQ n £ is independent of £ at 
fixed v for all modes analyzed here, where £ < 900 and vjL > 4.2 /xHz. Thus, the 
perturbation causing the frequency shift is probably localized in the region 4 Mm 
below the solar surface or less. Using the same reasoning, the thickness of the near- 
surface region where the uncertainties in the physics of the model are confined can 
be inferred. In this case, the set of modes where 5v n ,iQn,i is independent of £ at fixed 
v is such that v/L > 12 /iHz and £ < 370 using the high-£ frequency determination 
(Section 5.2). Including higher-degree modes, 5v n ^Q n i is not a function of frequency 
alone anymore (see, for example, Figure 3 in Rabello-Soares, Basu, Christensen- 
Dalsgaard, et al., 2000). As a result, the frequency uncertainties in the model are 
probably largely restricted to the layers 18 Mm below the surface. 



7. Conclusion 

In the determination of unbiased high-degree mode parameters, the instrumental 
characteristics must be taken into account in the image spatial decomposition or 
in the leakage matrix calculation itself to obtain a correct estimation of the relative 
amplitude of the spatial leaks. Among the instrumental characteristics analyzed here, 
the image scale is the one that affects the parameter determination the most. The 



p. 26 



Analysis of MDI High-Degree Mode Frequencies 



image scale is the ratio of the image dimensions observed on the CCD detector and 
the dimensions in the actual Sun. An error in the image scale introduces an error in 
the estimated central frequency which increases with the mode frequency. A 0.27% 
error in the image scale would shift the estimated central frequency by as much as 11 
/iHz at 5 mHz. The radial distortion also has an important effect which is expected 
since it is similar to an image scale error. An instrumental property not taken into 
account here (due to a lack of a good estimation) that could have an important effect 
on the measured parameters is an azimuthally varying PSF. 

The applied ridge-to-mode correction recovers frequencies at moderate degree that 
differ from the assumed corrected values by 1 /iHz or less depending on the mode 
frequency. The fitting uncertainty of the recovered frequencies is in the range 0.07- 
0.18 /iHz. The agreement for the ai, a2, and 03 coefficients is very good, except 
maybe for /- and pi-mode a\ coefficients, their mean difference with respect to the 
assumed correct values is, respectively, three and two normalized by the ridge fitting 
uncertainties. 

At high degree, the differences between our frequency determination and theo- 
retical frequencies for the p modes shows the same general variation with degree as 
the results obtained with ring analysis. For n < 5, the global and ring analysis agree 
within 6 /iHz. The high-degree /-mode frequencies obtained using ring analysis, 
like previous observations (Duvall, Kosovichev and Murawski, 1998 and references 
within), are substantially lower than the theoretical frequencies. Surprisingly, the 
/-mode high-£ set frequencies agree well with the model frequencies (within 3 /iHz) 
whereas the ring-analysis frequency differences can be as large as 13 /iHz for £ > 700 
modes. The implications of the high-/ frequencies and splitting coefficients on the 
solar structure and rotation will be addressed in a future paper. 

As noted by other authors for low- and moderate-degree modes {e.g., Libbrecht 
and Woodard, 1990), the frequency shift induced by the solar cycle scales well with 
the mode inertia. We extended this analysis to high-degree modes and found out that 
scaling with the mode inertia normalized by the inertia of a radial mode of the same 
frequency follows a simple power law (given by Equation 4) with one exponent at all 
frequency ranges, where the / and p modes are fitted independently of one another. 
The exponents obtained using four different solar indices agree within their fitting 
uncertainty, where: jf = 1.29 ±0.07 and j p = 3.60 ±0.01. The /-mode exponent is 
less than half of the p-mode value. The fundamental mode of solar oscillations has 
essentially the character of surface gravity waves and, contrary to the p modes, it is 
essentially incompressible and independent of the hydrostatic structure of the Sun. 
Hence, it is not a surprise that these different types of modes have different exponents. 
Due to their different properties, it is also very likely that different physical effects are 
responsible for their frequency variation. Accondingly to Dziembowski and Goode 
(2005), for the / modes, the dominant cause of frequency shift is the variation of 
the subphotospheric magnetic field. For the p modes, it is the decrease in the radial 
component of the turbulent velocity in the outer layers during the increase in solar 
activity, which is accompanied by a decrease in temperature (due to a decrease in 
the efficiency of convective transport). At low frequency (y < 2.3 mHz), the p-mode 
frequency shifts have a different behavior than at high-^: a step (with the same 
exponent 7p) or, as found by other authors, an exponent twice as large as the one 
at high-z/. Low-frequency modes have a large probability of being uncorrelated with 
the solar cycle, which was not taken into account in the case where a large exponent 
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was estimated. The /-mode frequency shifts also have a different behavior around 
1.7 mHz: they increase abruptly by an order of magnitude. 

Modes with frequency around 3 mHz have the smallest probability that their 
variation is linearly uncorrelated with the solar index, while modes with v < 2.5 
mHz or v > 4.5 mHz have the largest probability of being uncorrelated. A large 
probability (P u ) of being uncorrelated does not necessarily means that a given 
mode is not physically correlated with the solar cycle, instead it could be due to 
uncertainties in the measurements, a low signal-to-noise ratio. The logarithm of P u 
is well correlated with the logarithm of the relative uncertainty of the estimated 
frequency shift Sv e , the Pearson correlation coefficient is 0.71 for medium-£ modes 
and 0.54 for high-^ modes. If a given mode has a large probability of being linearly 
uncorrelated, it is expected that the linear fitting of its frequency shifts will have a 
large uncertainty, hence the high correlation coefficient between P u and the estimated 
frequency shift uncertainty. However, it raises the question of what could be the 
physical process that would make those modes less sensitive to solar activity. For 
a given £, the upper reflection point for lower-frequency modes is deeper in the 
Sun than for high-frequency modes. If the perturbation layer causing the frequency 
shift is above the upper turning point of the mode, it would not be affected by 
the solar cycle. Accordingly to model S of Christensen-Dalsgaard, the depth of the 
upper turning point increases sharply with decreasing frequency below 2.3 mHz. 
The upper turning point for a radial mode with v = 2 mHz is 0.5 Mm deeper 
in the Sun than a three-millihertz mode (from Figure 2 in Chaplin, Appourchaux, 
Elsworth, et al, 2001). At high-frequency, the observed frequency shift seems to 
suddendly drop to zero. However, there is no agreement on the exact frequency that 
this happens, the observed values range from 3.7 to 5 mHz. The frequency-shift 
falloff is explained by an increase of chromospheric temperature and magnetic field 
at solar maximum (Jain and Roberts, 1996 and references within). In the presence 
of an inclined magnetic field, high-frequency modes tunnel through the temperature 
minimum and are particularly sensitive to changes in the chromosphere, which are 
expected to be well correlated with solar activity. 
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